*Test the mean rank using the JTPA training data - Table 3, Panel A: Earnings, II: Mean Test
*Please make sure `ivqte' and `moremata' are installed
cd "C:\Users\XXX"
*Replace `Gender' by male or female
log using ".\Dropbox\Rank Invariance\Empirical application\JTPA\\Results\Table3_IIA_Gender.smcl", replace
clear all
set more off
set matsize 5000
version 12

	capture program drop test_mean
	program define test_mean, eclass 
	tempname diff qte qt_c qt_t
	matrix `diff'=J(1,numl,.)
	tempvar y2 rank
	qui gen `y2'=y
	qui replace `y2'=0 if t==1
	qui ivqte y (t=z), quantiles(0.01(0.01)0.99) 
	matrix `qte'=e(b)
	qui ivqte `y2' (t=z), quantiles(0.01(0.01)0.99)    // use y2 to trick ivqte estimate quantile under control
	matrix `qt_c'=-e(b) // notice the negative sign
	matrix `qt_t'=-e(b)+`qte'
	forv q = 1(1)99 {
	qui gen uind`q'=y>=`qt_c'[1,`q'] if t==0
	qui replace uind`q'=y>= `qt_t'[1,`q'] if t==1
	}
	egen `rank'=rowmean(uind1-uind99)	
	foreach l of global level {
			local lind=`lind'+1
			qui sum `rank' if z==0 & x==`l'
			scalar temp=r(mean)
			qui sum `rank' if z==1 & x==`l'
			matrix `diff'[1,`lind']=r(mean)-temp
	}
	cap drop uind*
	ereturn post `diff'
	end 

use ".\Dropbox\Rank Invariance\Empirical application\JTPA\JTPA_ContinousAge.dta", clear
gen y=earn
gen z=assignmt
gen t=training
qui sum y
local obs=r(N)
replace hsorged=round(hsorged)
replace married=round(married)
replace wkless13=round(wkless13)
*Column I: Female
local male=0
**Column I: Male
*local male=1

*cov =0/1 controls for whether using covariates in the first-stage QTE estimation.
local cov=1
local mincell=5
local bsreps=2000

*age is in 6 categories, which are represented by 5 age dummies
local cov1 "hsorged  black   hispanic  married  wkless13  age2225 age2629   age3035  age3644   age4554 "
local cov0 "hsorged   black    hispanic   married   wkless13 afdc age2225    age2629    age3035 age3644   age4554 "
qui tostring `cov1', replace
qui tostring `cov0', replace
*sex==1 for male, sex==0 for female	
gen x=hsorged+black+hispanic+married+wkless13+age2225+age2629+age3035+age3644+age4554 if sex==1
replace x=hsorged+black+hispanic+married+wkless13+afdc+age2225+age2629+age3035+age3644+age4554 if sex==0

qui destring x, replace	
qui destring `cov1', replace
qui destring `cov0', replace

if `male'==1&`cov'==1{
keep if sex==1
global cov="`cov1'"
}
else if `male'==0&`cov'==1 {
keep if sex==0
global cov="`cov0'"
}
else if `male'==1&`cov'==0 {
keep if sex==1
global cov=" "
}
else {
keep if sex==0
global cov=" "
}

qui tab x, sort
qui levelsof x, local(level)
local numl0 : word count `level'
scalar numl0=`numl0'
scalar di numl0

gen flag=1
foreach l of local level {
	qui count if (x==`l' & z==0)
	scalar nullcell0=r(N)
	qui count if (x==`l' & z==1)
	scalar nullcell1=r(N)
		if (nullcell0<`mincell'|nullcell1<`mincell') {
		qui replace flag=2 if x==`l'
		}
local minl=`l'
}

qui sum flag if flag==2
	if r(N)==0 {
	*If no cells are dropped at all, need to drop one out of the total number of unique values of X
	drop if x==`minl' 
	}

qui levelsof x if flag==1, local(level)

global level= "`level'"	
local numl: word count $level
scalar numl=`numl' 
*numl is the actual number of X values used
scalar di numl

	set seed 123
	timer on 1
	if `bsreps'==0 {
	di in red "Need to set the number of bs replications."
	exit
	} 
	else {	
	qui test_mean
	matrix dm=e(b)
	bootstrap _b, reps(`bsreps'): test_mean
	matrix hatv1=e(V)
	scalar df1=numl-diag0cnt(hatv1)
	matrix wald1=dm*invsym(hatv1)*dm'
	scalar bsstat1=wald1[1,1]
	scalar Pval1=1-chi2(df1, bsstat1)
	scalar dec1=(wald1[1,1]>invchi2(df1,0.95))
	}  	
scalar mincell=`mincell'

scalar di mincell
scalar di numl0
scalar di numl
scalar di df1
scalar di bsstat1
scalar di Pval1
scalar di dec1

timer off 1
timer list 1
scalar time=r(t1)/60
scalar di time


log close

	


